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Abstract 

Using methods from effective field theory, we develop a novel, systematic framework for 
the calculation of the cross sections for electroweak gauge-boson production at small 
and very small transverse momentum qt, in which large logarithms of the scale ratio 
Mv/qt are resummed to all orders. These cross sections receive logarithmically en- 
hanced corrections from two sources: the running of the hard matching coefficient and 
the collinear factorization anomaly. The anomaly leads to the dynamical generation 
of a non-perturbative scale q^: ~ My e~'^°°^*/"°*-^^^\ which protects the processes from 
receiving large long-distance hadronic contributions. Expanding the cross sections in 
either as or qt generates strongly divergent series, which must be resummed. As a by- 
product, we obtain an explicit non-perturbative expression for the intercept of the cross 
sections at qt = 0, including the normalization and first-order as{q*) correction. We 
perform a detailed numerical comparison of our predictions with the available data on 
the transverse- momentum distribution in Z-boson production at the Tevatron and LHC. 



1 Introduction 



In collider processes with several disparate scales, fixed-order perturbative expansions in QCD 
become unreliable since higher-order corrections are enhanced by large double logarithms of 
scale ratios. The classic example of such a multi-scale process is the Drell-Yan production 
of electroweak gauge bosons with transverse momentum qt much smaller than their mass. 
The leading logarithmically-enhanced corrections in this kinematic region were resummed in 
[1, 2, 3], and an all-order formula for the resummed cross section at small qr was obtained in 
the seminal work [4]. The region of small qt is of great phenomenological importance, since it 
has the largest cross section and is used e.g. to extract the W^-boson mass and width. In the 
related process of Higgs-boson production via gluon fusion, the region of small qx is important 
because one usually vetoes hard jets in order to enhance the signal over background ratio. The 
traditional resummation approach of [4] suffers from the presence of singularities arising from 
integrals over the Landau pole of the running coupling constant, and hence a prescription 
is required to regularize the integral. In practical applications, the integration is cut off 
at large xt values, and to account for the missing contributions a non-perturbative model 
function is employed. For transverse momenta in the perturbative domain these long-distance 
contributions are formally power suppressed, but it is irritating that an explicit prescription 
for how to deal with them is needed even for qt values deep in the perturbative regime. The 
explicit cut-off also makes it difficult to perform the matching onto fixed-order computations, 
since cut-off effects persist even when the formula is evaluated at large qt ~ My, where My 
denotes the mass of the Drell-Yan object. 

In a recent paper, we have derived a novel form of the all-order factorization theorem for 
Drell-Yan production at small transverse momentum, Aqcd <^ Q't ^ My [5], using methods 
from effective field theory. The cross section is expressed as a product of a g^-dependent hard 
function with a convolution of two transverse-position (xt) dependent parton distribution 
functions (PDFs). The form of the factorization theorem is affected by an anomaly of the 
effective Lagrangian, which implies that the naive factorization of two coUinear sectors valid 
at the level of the classical Lagrangian is broken by quantum effects. This gives rise to 
an anomalous, g^-dependent factor under the convolution integral in space. Consistency 
conditions ensure that this factor is a pure power of g^, with an XT^-dependent exponent. In 
the short-distance region, for xt <^ ^qcd' product of transverse-position dependent PDFs 
can be expanded in terms of standard PDFs convoluted with perturbatively calculable kernel 
functions. The resulting effective field-theory formula for the resummed cross section is free 
of Landau-pole singularities and therefore per se does not require non-perturbative modeling. 

Naively, assuming that transverse momentum and transverse displacement xj- are conju- 
gate variables satisfying QtXt ~ 1, the factorization formula derived in [4] and [5] only applies 
as long as qr ^ Aqcd is in the perturbative domain. However, as early as in 1979, Parisi 
and Petronzio have argued that for asymptotically large momentum transfer, the Drell-Yan 
cross section at very small and even vanishing transverse momentum of the lepton pair can 
be calculated in resummed perturbation theory [2]. Using methods from effective field theory, 
we develop a novel framework for the systematic short-distance calculation of the Drell-Yan 
cross section at small and very small qr- We show that if the mass My of the Drell-Yan 
object is so large that the scale g* ~ My exp ~ (^^Cr+PoW (Mv) ™ perturbative domain 
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(here R = F,A denotes the color representation of the partons creating the object), the cross 
section can be calculated using renormalization-group (RG) improved perturbation theory for 
arbitrarily small values of qt, up to power corrections controlled by Aqcb/q*- The scale 
emerges dynamically, and it screens the cross section from receiving long-distance contribu- 
tions at leading power. We distinguish three regions of transverse momentum, which require 
different classes of terms to be resummed. In the region of large qx, approximately q^ > 15 GeV 
for the case of electroweak gauge-boson production, a fixed-order perturbative calculation is 
justified. In an intermediate region, the perturbative series exhibits large logarithms as well 
as a strong factorial growth of a certain set of higher-order corrections, and both must be 
resummed in order to obtain reliable results. Finally, in the region of very small < g^,, 
the power counting of the perturbative expansion must be modified and certain higher-order 
terms must be resummed at all orders, even though they are naively suppressed by powers of 
as- We set up our formalism in such a way that a single formula automatically interpolates 
between the three regions. We provide all ingredients to perform the resummation of large 
logarithms at next-to-next-to-leading logarithmic (NNLL) order. In the low-gr region this 
requires retaining terms of up to four-loop order in some coefficient functions. These terms 
can be related to known anomalous dimensions using RG equations. 

As a by-product of our analysis, we derive for the first time an explicit formula for the inter- 
cept of the differential cross section da/dq^ at vanishing transverse momentum, in which the 
normalization and first-order perturbative correction in as{q*) are included. We show that in 
the region of very small q^ a usual operator-product expansion (OPE) of the cross section is in- 
applicable due to a very strong asymptotic divergence of the twist expansion. The dependence 
on qT and the sensitivity to long-distance hadronic effects can only be assessed by resumming 
the OPE. For the case of the qr dependence such a resummation is implicit in our approach. 
Long-distance corrections, which could affect the matching of transverse-position dependent 
PDFs onto standard PDFs, can at present only be modeled using a phenomenological function 
fhadiixT-^Np)- Our formalism provides a convenient framework for implementing such a form 
factor. Using different forms of model functions, we find that to a good accuracy the shape 
of the cross section is only affected by the coefficient of the first term in the expansion of the 
form factor in transverse separation, fhadiixr-^Np) = 1 — A^p + • • • • The main effect of 
long-distance corrections is to shift the peak of the da/dqr distribution for Z bosons produced 
at the Tevatron and LHC by roughly about Anp. We do not confirm the rather pessimistic 
statements about the vafidity of the short-distance analysis of the transverse-momentum dis- 
tribution in the region of very small qr made in [4] , where it was concluded that a Drell- Yan 
mass as large as 10^ GeV would be required to keep power corrections to the intercept at 
an acceptable level (below 20%). Instead, we find that the theory developed in this paper 
works well for W and Z production at hadron colliders, and it will work even better for the 
production of Higgs bosons via gluon fusion or the production of new heavy particles such as 
W or Z' bosons or slepton pairs. 

In Section 2 we briefly review results from our previous work [5] and reorganize them in a 
way suited for our discussion. In Section 3 we then explain the systematics of the resummation 
scheme valid at very small g^, where the emergent non-perturbative short-distance scale 
protects the cross section from receiving long-distance contributions at leading power. The 
modified power counting required in this region is introduced and used to derive the relevant 
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expansions of the various coefficient functions in the factorization formula for the differential 
cross section. We also explain the relevance of the three regions in transverse-momentum space 
mentioned above. Section 4 is devoted to a study of the qt distribution near the origin. We 
derive an explicit formula for the cross section da/dq^ qt — including the normalization 
and ffist-order correction in as{q^). We also explain why the OPE breaks down for very small 
qt- Detailed systematics studies of of our results arc performed in Section 5, while Section 6 
is devoted to comparisons with Tevatron and LHC data on Z-boson production. We conclude 
in Section 7. Technical details of our calculations are described in four appendices. 



2 Factorization and resummation 

In a recent paper [5] , we have analyzed the Drell-Yan process using methods of effective field 
theory. In the kincmatical region where Aqcd ^ (/t ^ Mz, the double differential cross 
section for Z-boson production was shown to obey the factorization formula 
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where g1 ^ denote the Z-boson couplings to quarks qL,R in units of e, (f)i/N{z,fi) are standard 
PDFs, and ^1,2 = where r = (M| + qj')/s and y is the rapidity of the Z boson in 

the laboratory frame. A corresponding formula for the differential cross section da/dq^, along 
with explicit expressions for the weak charges, is presented in Appendix A. The function Cy 
is a short-distance ( "hard" ) Wilson coefficient arising in the matching of the weak currents of 
the Z boson onto the leading-power current operator in soft-collinear effective theory (SCET) 
[6, 7] . The perturbative ( "coUinear" ) kernel functions Cgg^ij are given by (we denote q^ = —q\ 
and = —x\) 



2 A/f2 \ -Fqq(L±,as) 

(2) 

X Iq<-i{zi, L±, Us) Iq^j{z2, L±, Ug) , 



where 



"^^^^ ^ -ln%!, 6o = 2e--. (3) 



47r ' - 6g 



The functions li^j arise in the matching of transverse-position dependent PDFs Sj/jv defined 
in [5] onto ordinary PDFs, 

Bi/N{^, xl,iJ,) = ^ — Ii^j{z, xl, 11) (t)j/N{i/z, n) + 0(A^cD ^t) > (4) 
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which is vahd at small transverse separation x"^ ^ ^qcd- The factorized cross section (1) 
receives power corrections in the two small ratios q^/M^ and Aqcd/Qt^ where the latter ones 
enter via (4). These power corrections will not be indicated explicitly in most of our equations. 

In the factorization formulas (1) and (2), the dependence on the two disparate scales Mz 
and Qt is factorized explicitly in Xt space. Note the unusual fact that dependence on the 
high scale Mz enters in two places: via the hard matching coefficient Cy, but also via an 
XT-dependent power of Mz under the Fourier integral in (2). The latter effect is due to the 
coUinear factorization anomaly discovered in [5] (see also [8]), which results from the fact 
that the naive factorization property of the classical effective Lagrangian of SCET (i.e., the 
property that different sectors of the effective theory do not interact with one another) is in 
some cases spoiled by quantum effects. This anomaly introduces a power-law dependence on 
Mz in the kernel functions in xt space. The anomalous exponent Fqq is constrained by the 
non-abelian exponentiation theorem [9, 10]. As long as x"^ <C Aqq^, it can be calculated in 
perturbation theory, and it is presently known to two-loop order [5]. A simpler example for 
the occurrence of the coUinear anomaly is the Sudakov form factor of a massive vector boson, 
for which the anomalous power-law dependence on the gauge-boson mass was derived in [11]. 
More recently, it was shown that also the jet-broadening distribution in e~^e~ annihilations is 
affected by a factorization anomaly [12]. 

The factorization formulas (1) and (2) have been shown in [5] to be equivalent, to all orders 
in perturbation theory, to an expression for the resummcd cross section derived in a seminal 
work by Collins, Soper, and Sterman (CSS) [4]. An interesting effect of the anomalous terms 
is that they give rise to additional contributions to the resummation exponents A and B in 
the CSS approach, which therefore do not agree with the corresponding quantities for soft- 
gluon resummation. In particular, starting at three-loop order the coefficient A is longer equal 
to the cusp anomalous dimension F^^p. We note at this point that recently an alternative 
SCET-based resummation scheme for Drell-Yan production at small Qt was presented [13, 14], 
in which only the large logarithms contained in the hard function |Ci/p in (1) are resummed. 
This approach does not address the resummation of the large logarithms steming from the 
coUinear anomaly, which reside in the kernel functions Cqq^ij in (2). The method of [13, 14] 
therefore does not provide a consistent resummation scheme already at NLL order. While 
these authors correctly reproduce the next-to-leading logarithms in the cross section itself, 
this is not sufficient. For observables subject to Sudakov double logarithms, it is essential 
that the counting of logarithms is performed in the exponent and not at the level of the cross 
section. If the NLL terms are not exponentiated, higher-order terms in perturbation theory 
become arbitrarily large in the region where the logarithms are of 0{l/as)- In Appendix D 
we discuss the systcmatics of the expansion in more detail and give explicit expressions for 
the fixed-order expansion to 0{a^g) of our result for the cross section. We also show explicitly 
which of these terms contain the large logarithms not accounted for by the resummation of 
the hard function. 

The resummation of large logarithms in the factorization formula (1) is accomplished by 
evolving the hard matching coefficient Cy to a scale determined by the average transverse 
separation, fi ~ {x^^) <^ M^, in a sense described in more detail in the next section. With 
this scale choice the logarithm Lj^ in (3) is of C(l), and hence the functions F^^ and Iq^i can 
be calculated using fixed-order perturbation theory. All large logarithms are then contained 
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in the coefficient Cy(—M'^,fi) and in the anomalous factor (x'^M'^)^^''^ . The RG evolution 
equation of the hard matching coefficient at time-like momentum transfer is of the Sudakov 
type and reads [7] 

^2 



d 



Cvi-q\l^) 



Cv{-q\l^), (5) 



dln/j, 

where T^gp is the cusp anomalous dimension in the fundamental representation and 7^ denotes 
the anomalous dimension of a coUinear quark field in SCET. These quantities are known to 
three-loop order. The explicit form of the solution to this equation up to next-to-next-to- 
Icading order (NNLO) in RG-improved perturbation theory has been discussed in detail in 
[15]. The advantages of using a time-like scale choice {fi\ < 0) for time-like processes such 
as Drell-Yan production were emphasized in [16, 17, 18]. In Appendix B, we compile the 
NLO expression for the hard matching coefficient needed for our analysis. RG invariance of 
the cross section (1) requires that the cusp logarithm 2r^gplnM| resulting from the scale 
variation of the hard function be compensated by a corresponding term in the scale variation 
of the kernels Cgg^^ij. This is ensured by the RG equation [5] 

Fgg{L^, as) = 2rf„gp(a,) (6) 

for the anomalous exponent. Even though it is not required for the resummation procedure, it 
will be important to also consider the evolution equations for the kernel functions Iq^-i- They 
are given by 

Ig^i{z, a,) = [rf„gp(as) Ls_ - 27^(05)] Ig^i{z, a^) 



/" du 

- / — Ig^j{u,L±,as)Vj^i{z/u,as) , 



where Vj^^i are the usual DGLAP splitting functions. The first term on the right-hand side 
implies that the functions Ig^i exhibit double logarithmic dependence on in the exponent. 
It will be important for our purposes to factor out these terms, and this can be accomplished 
by rewriting 

IgU^, L^, as) = e'*^^^^''^^) Lx, a,) , (8) 

where ^ 

hpiL^, as) = rf„sp(as) Lj^ - 2-f'^{as) . (9) 

We choose to define hp{0,as) = 0, so that hp{Lj_,as) contains logarithms of Lj_ only. The 
new functions Ig^i now evolve exactly like the usual PDFs (but with the opposite sign in front 
of the DGLAP sphtting functions), while hp contains all double- logarithmic terms. We can 
now rewrite the hard-scattering kernels from (2) in the form 



Cgg^ij{zi, Z2, Qt, Ml, jj,) ^ - J dxT XT MxtQt) exp griM^, ji, L±, as) 

X Iq-<^i{Zl, L±, as) Iq^j{z2, L±, Og) , 



(10) 
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where 



9F{Ml,iJ,,L_L,as 



+ L_L ) Fgg{L±, Qs) + 2hF{L±, Qs 



(11) 



As mentioned above, our strategy will be to choose the factorization scale fx ~ {xj^^} <^ Mz, 
in such a way that the functions Fgg and hp, as well as the kernels Ig^i, can be evaluated 
in fixed-order perturbation theory. When using the resulting expressions to write down the 
perturbative expansion of the exponent g'F, we treat In (M|//i^) as a large logarithm and count 



7] = TjF^Mz.ij) = Tq a. In 



TT 



(12) 



as an variable. Here = ACp is the one-loop coefficient of the cusp anomalous dimen- 
sion in the fundamental representation of SU(A^c)- We then obtain 



9F{r), L±, as) = -r]L± - as 



(r^ + ^^o) -f + (27o' + VK) + rid2 



+ 0{al) , 



(13) 



where 7o = —SCp and /5o = y Ca — 3 TpUf are the one-loop coefficients of the quark anoma- 
lous dimension and /3-function, and with a slight abuse of notation we have eliminated the 
arguments M| and in gi^ in favor of rj. The quantities 



K 



rf 



67 20 
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contain some two-loop information. These particular ratios are the same for any representation 
of the gauge group. Finally, the one- loop expressions for the kernel functions Ig^^i in (10) read 



Iq^i(z, L_L, as) ^ 5(1- z) 5gi - as 



+ 0{al) , 



where 



riU-) = , PW,(.) = 4Tp [z' + (1 - z)'] 



are the one-loop DGLAP splitting functions, and the remainder functions 



(15) 



(16) 



TT 



2{l-z)--S{l-z) 



can be extracted from results obtained in [5]. 



TZg^giz) ^ ATf z{l - z) 



(17) 



3 Scale setting and systematics of the expansion 

A reliable evaluation of the factorization formula (1) requires that the factorization scale fi be 
chosen such that the logarithm L± = ln(a;fn/i.^/6o) entering the various coefficient functions in 
(10) is a small quantity. The most naive choice would be to set ~ ^q/^^t inside the integral 
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over the Bessel function, so that L± is a small logarithm for any choice of x^- There are 
several disadvantages to that treatment. First, since xt is integrated over all possible values, 
there is no clear meaning to the scale /i in the sense of a physical, characteristic scale of the 
process. Second, setting the scale inside the integral imphes that the integration unavoidably 
hits the Landau pole of the running coupling, giving rise to ambiguities in the numerical 
results. In the spirit of effective field theory, the scale /i should correspond to a physical scale 
in the underlying factorization theorem. In our case, where we have chosen to evolve the 
hard function down to a typical coUinear scale, the requirement is that /i should be set such 
that the perturbative series for the kernels Cqq^^ij has a well-behaved perturbative expansion. 
Obviously, this requires that on average the XT-dependent logarithm L±_ is small, so that the 
perturbative expansions of the functions gp and Jg^j are well behaved. 

The factorization formula (1) has been derived under the assumption Aqcd ^ Q't ^ Mz- 
It neglects power corrections of order g|-/M|, and also higher-order terms ~ h?^QYyx\ in (4). 
The first type of power corrections can be included by matching our formulas to the known 
fixed-order results for the differential cross section (see below). Naively, one would expect that 
the transverse momentum and transverse separation x^ are conjugate variables satisfying 
qrpXT ~ 1, in which case the second class of power corrections would scale like ^^cTi/qr- 
While this is sometimes true, the general situation turns out to be more complicated. After 
integration over the factorized dependence on Mz and qx in (2) gets intertwined in a 
complicated way, and this gives rise to the peculiar effect that the two scales qx and xt 
decouple for very small qT as long as Mz is sufficiently large. This explains an observation 
made long ago by Parisi and Pctronzio, who found that the intercept of the differential Drell- 
Yan cross section da/dq^ at gr = can be calculated using short-distance methods provided 
that the mass My of the Drell-Yan pair is asymptotically large [2]. Note that without the 
coUinear anomaly (i.e., for Fqq — 0) such an effect could not arise, since then the scales and 
Xt would trivially be related by gj_ ■ a;_|_ ~ 1 due to the Fourier integral in (2). The interesting 
interplay beteween the scales My and q^p is less obvious in the traditional formalism, which 
implicitly adopts the choice ~ bo/xr- This choice eliminates the anomaly at NLL, but the 
My dependence then enters via the hard function, which for /i bo / xt is a, non-trivial part 
of the Fourier integral. 

To analyze the Fourier integral in the effective theory, consider the leading-order approxi- 
mation for the kernels Cqq^ij in (10), where we set — )■ and only keep the first term in the 
perturbative series for the exponent gp in (13). This gives rise to the integral 

i/V../o(..,.)e-^^(^)"^. (18) 

The integral over the Bessel function converges at large distances {xt oo) only ii r] > |, 
but the integral can be defined by analytic continuation for all values < r] < 1, and it 
is equal to the original Fourier integral (2) in the distribution sense. Provided rj is in this 
range, the result confirms the expected scaling xt ~ Qt^ (since e"^^-^ oc (x|./x^)~^) modulo a 
numerical (?7-dependent) factor. Contributions from large values xt ^ q^^ are suppressed due 
to the rapid oscillations of the Bessel function, while contributions from small values Xt ^ Qt^ 
are phase-space suppressed. Additional powers of L± under the integral in (18), which enter 
when the 0{as) corrections in (10) are included, can be generated by taking derivatives with 
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respect to rj. It follows that the scale choice ji ^ indeed eliminates large logarithms in the 
(Fourier-transformed) coefficient functions, as might have been expected from the beginning. 

The situation changes, however, at small transverse momentum. There are two possibilities 
that must be differentiated. Consider first the case where 77 defined in (12) remains smaller 
than 1 when fi approaches the non-perturbative domain. In this case the above discussion 
remains valid, but the perturbative calculation of the kernels Ig^i is no longer under theoretical 
control. Indeed, in this case one must resort to a formula for the cross section in terms of 
the transverse-position dependent PDFs Bi/N{^,x'^, fj.) in (4), which are then genuine, non- 
perturbative objects. But there is also a second possibility, that 77 reaches 1 for values of /i 
that are still in the perturbative domain. We denote by the value of /i where this happens, 
i.e. 

= "^'"^ " (" (rf + ftK(M.) ) - 

where in the last step we have used the one-loop approximation for the running coupling. 
Solving the first equation numerically, one finds ~ 1.88 GeV in the present case, which is 
indeed a reasonable short-distance scale. For 1] > 1, the integral in (18) is ultra-violet (UV) 
divergent for Xt 0, since the integrand then approaches (pl/ i^PY xlf"^^ . It is then necessary 
to keep higher-order terms in the perturbative series for the exponent gp in (13). Indeed, 
since the quadratic term in L\ has a negative coefficient, it provides a gaussian weight to the 
integral which cuts off the divergence at xt — > 0. But at the same time, the quadratic term 
also provides a regulator for the infrared (IR) region of very large xt, which is in addition 
to the oscillating behavior of the Bessel function. For small q^ < q*, this gaussian fall-off 
is the dominating factor, which prevents that xt can become arbitrarily large. Remarkably, 
this implies that (xt) decouples from q^^ and stays in the short-distance domain even in the 
extreme case where qr is taken to 0. Using (13) and changing variables from xt to i — L±, 
we obtain 



1)2 

^?5-^jj L„-^n = / d,£ exp 



•2 

2\ 



iz2,i,as). (20) 



The integrand features a gaussian peak at 



-peak 



(21) 



with a width proportional to 1/ ^/a^. The condition that at the peak the logarithm i — L± 
should be an 0{1) quantity implies that 1 — rj — 0{as), implying that the factorization scale 
must be chosen in the vicinity of the scale q^ in (19). In other words, once the scale fj, reaches 
g*, the value of jj, that keeps the logarithms L± small decouples from qT and stays near g* even 
when qT ^ 0: 

fi ~ (xy^) ~ max{qT, q*) ■ (22) 

In our numerical work, we will use ^ = qr + Q* as the default choice. We recall that the emer- 
gence of the scale q^ in (2), below which the scaling of x^^ decouples from qx, is a consequence 
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of the collinear anomaly, which is responsible for the anomalous dependence on Mz in (2). 
Provided the mass of the Drell-Yan boson is large enough that 3> Aqcd, the transverse- 
momentum distribution is protected from long-distance physics even for arbitrarily small qt- 
The resummed perturbative series for the cross section generates the scale 5* dynamically, and 
even though this is a short-distance scale, it is related to the boson mass Mz in a genuinely 
non-perturbative way. 

The above discussion shows that we must distinguish two regions of transverse momenta. 
For qx 3> the Bessel function regularizes the UV region, and the scale choice /i ^ qr prevents 
that the logarithms L± give rise to large perturbative corrections. It is then consistent to count 
these logarithms as L± ~ 1 and construct the perturbative series as a series in powers of a^. For 
qr ^ Q* the situation is different. Even though the scale choice ^ ^ q^ ensures that = (9(1) 
at the peak of the integrand, the gaussian weight factor allows for significant contributions to 
the integral over a range of Lj_ values with width proportional to see (20). It follows 

that for very small gr the resummation procedure must be reorganized, using the modified 
power counting L±_ ~ l/^/oJ. This implies that single-logarithmic terms (a^Lj^)"^ ~ al^"^ are 
always suppressed, whereas double-logarithmic terms (agL^^)"" ~ 1 are unsuppressed and must 
be resummed to all orders. To keep track of this fact, we introduce an auxiliary expansion 
parameter e (which at the end is set to 1) and assign the power counting ~ e and L±_ ~ e~^/^. 
In Appendix C, we use the recursive solutions to the RG equations (6) and (9) to determine 
all terms in Fqq and hp that contribute up to 0{e) to the exponent gp defined in (11). This 
involves some four- loop contributions to Fqq and some three-loop contributions to hp, which 
however can all be expressed in terms of one- and two-loop coefficients of the anomalous 
dimensions and ^-function. The resulting expression is 

as (r^ + 77/3o) ^ ' 



2 







a, (27o^ + rjK) + a] (F^ + 77^0) ^0 3 



,1/2 



(23) 



as vd2 + al (XF^ + 270^^0 + ^7 (^1 + 2X/3o)) ^ + (F^ + rj^o) Pi - 

- e»(e^/2) , 

which is more complicated than the naive perturbative expansion in (13). The auxiliary 
parameter e counts the order in resulting (for qT <S q*) after the xt integral in (10) has 
been performed. The two terms given in the first line are unsuppressed and must be kept in 
the exponent of the integral in (10), whereas the remaining terms can be expanded in powers 
of e^/^. The resulting integrals over the Bessel function in (10) can readily be evaluated 
numerically. An efficient way of doing this is to use that Jo^XTqr) = f Ivn KQ^—ixTqr) and to 
perform a contour rotation from xt — > ixt- 

It is interesting to ask whether the "protective behavior" in both the UV and IR regions 
provided by the gaussian terms in the exponent could be upset at yet higher orders in the 
perturbative expansion of the exponent. It is not difficult to show that the highest-order 
logarithmic terms in (23) are given to all orders by — a"+^ (F^ -|- r]j3o) /3q (n -|- 2) with 
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n > 0. Provided that the series is truncated at an integer order in e, in which case the largest 
n value is even, these terms are negative, and they guarantee that the integral over xt in (10) 
converges both in the UV and IR regions. When the expansion is performed in the exponent, 
the series of these highest logarithms does not exhibit a factorial divergence. On the contrary, 
the sum is convergent as long as L± < l/(/3oas)i which is always the case with our power 
counting L± ~ 1/ ^/aj- Hence, the modified power counting displayed in (23) indeed provides 
a consistent, well-behaved expansion scheme for the integral over the Bessel function. 

We finally give the expressions for the coUinear kernel function arising with our modified 
power counting. We find 



Iq^i{z, L_L, as) = 5{1 - Z) 5qi 

asU^U^) + al [V,U^) - 2^o7^i^.(^)) ^ 



+ 



(24) 



which may be compared with (15). Here 



J z " 

3=(l,9 



(25) 



involve the convolutions of two DGLAP splitting functions. The resulting expressions are 
given in Appendix C. 

Formulas (23) and (24) are our main results. With the help of these expressions, large 

logarithms can be resummed at NNLL order all the way down to zero transverse momentum 
(always assuming that is in the perturbative domain, as is indeed the case for Z-boson 
production). For larger qx values the additional terms contained in (23) and (24) compared 
with (13) and (15) reduce to higher-order terms proportional to and a^, which can be 
neglected to the order we are working. Hence, our formula provides a smooth interpolation 
between the regions of small and very small q^. 

Interestingly, the additional terms needed at very low q^ are also important at qx > q*, as 
we will now show. We have argued above that the quadratic terms in L± need to be kept in 
the exponent as rj approaches 1, in order to regularize the integral over xt in the UV region. 
As long as r) is less than 1 (i.e. for q > q*), it would seem justified to expand out the 0{as) 
correction to the exponent gp in (13) in a power series in a^. Surprisingly, it turns out that 
this would be a bad idea. Closer inspection of the integral in (18) shows that even for (1 — 77) 
small but not zero, taking n derivatives with respect to t] (corresponding to an extra factor 
{—L±)"' under the integral) generates a contribution exhibiting factorial growth in n. For 
the special case /i — qr the series has been analyzed in [5], where it was shown that (setting 
a = 2as(ro -I- rj/3o) for brevity) 



1 



dxT Xt Jo^XTqr) e 



5:|^(-i)' 

n=0 



expanded 
A'=?T, r)<\ 



1 



-27B 



(1-77) 



2n+l 



+ A;„a" + C»(l-7y) L 



(26) 
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where the coefficients A;„ do not exhibit the strong factorial growth of the terms in the ffist 
sum. (Numerically, we find ko ^ 0.315, h ^ 0.243, k2 ^ -0.075, ^ -0.051, ki ^ 0.045, 

) While the series in (26) is badly divergent, the fact that it has alternating sign implies 

that it can be Borel summed, yielding 



1 r 

2 Jo 



a ^2 Borel sum 



e-2T^^ /TT 



- e 



1-Erf 



1—7] 



1-Erf 



(27) 



oo 



+ ^ V A;„a" + e»(l-77), 



n=0 

where Erf(x) is the error function. This expression provides an explicit (albeit approximate) 
result for the complicated integral on the left-hand side. Note that the spurious singularity 
at ?7 = 1 in the sum (26) has disappeared after Borel summation, but that the final answer 
depends in a highly non-perturbative way on the coupling constant ~ a. Even though a 
perturbative expansion in powers of a ~ can be constructed, its radius of convergence is 
zero. The lesson from this exercise is that even for values of 77 large but not very close to 1 it 
is obligatory to keep the 0{af,L\) terms in the exponent of the xt integral, as suggested by 
our modified power counting in (23). Expanding these terms out is only justified for rather 
large g^, for which r] = 0{as) is a small quantity. 

How can this strange behavior be understood? The integral over the qt spectrum from up 
to some value g™^ » Aqcd should be calculable in fixed-order perturbation theory, as is the 
integral over all possible values. Why do problems arise for smaller values of qt, which naively 
are still in the perturbative domain? As we will show in the following section, the intercept 
of the spectrum da/dq^ (or d'^ a / dq^dy) can be calculated using short-distance methods, but 
it exhibits an essential singularity at ttg = and hence does not admit any perturbative 
expansion. Since the spectrum is at least approximately given by a smoothly falling function, 
the area under it (the total cross section) is roughly proportional to the height (intercept) 
times a typical width. The requirement that the area be perturbatively calculable requires 
that not only the height, but also the characteristic width of the spectrum are governed by 
non-perturbative short-distance physics. 

In order to make these comments more precise and see the relevance of the different con- 
tributions discussed in this section, we show in Figure 1 the differential cross sections da/dq^ 
(left) and da/dqT (right) for Z-boson production (with subsequent decay to a lepton pair 
Z — 7- ^'^^^) at the Tevatron, using three different approximation schemes. The green bands 
show the cross sections at NNLL order in our improved expansion scheme based on the modi- 
fied power counting given by the e expansion discussed above ( "improved resummation" ) . The 
red bands show results obtained with a conventional power counting in a^, keeping however 
the quadratic terms of order asL\ in the exponent ("conventional resummation"). Finally, 
the orange bands correspond to the resummation scheme where all terms in the exponent 
are expanded out ("naive resummation"). Considering first the plot on the left, we see that 
as expected the improved resummation becomes important for gr ^ ~ 1.88 GeV, where it 
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qr [GeV] qj [GeV] 

Figure 1: Comparison of different expansion schemes. The bands result from varying fi by 
a factor two around the default value fi = qt + q*- The thick lines in the left plot show 
the default predictions. The wide orange bands are obtained using the naive resummation 
formula, which suffers from a factorial divergence. The red bands result if the divergent terms 
are resummed (conventional resummation). The dark green, densely hatched bands arise if 
one further resums the terms which are enhanced at very small qt (improved expansion). 

gives results that are lower and have a smaller scale uncertainty compared with the conven- 
tional resummation scheme. In contrast, the naive expansion, which is affected by the factorial 
divergence, leads to very large scale dependence, such that the cross section becomes negative 
if the scale is lowered by a factor of two from its default value, and even the default prediction 
shown by the solid thick orange curve becomes negative for qj- < 1.5 GeV. In the region of 
very small transverse momentum, the conventional resummation scheme is only valid to NLL 
accuracy, since it misses the additional higher-order terms which are enhanced in this region. 
As a consequence, the scale uncertainty near qt = is much larger in this case than for the 
improved resummation scheme. Interestingly, the effects of including the additional terms is 
also non-negligible at higher values of qt- While the results obtained in the conventional and 
improved resummation schemes are compatible within scale uncertainties, the improved result 
has a significantly smaller uncertainty, and the peak of the distribution is shifted slightly to 
the right. The reason that the effects of the additional resummation are still visible away from 
the end-point is the strong fall-off of the spectrum towards larger qt- 

In addition to the logarithmic terms which are resummed by our result, the cross section 
also contains regular terms, which can be obtained from a fixed-order computation of the 
spectrum. In order to capture both corrections, we combine our result with the fixed-order 
result for the qr spectrum. To avoid double counting of the logarithmic terms, we need to 
subtract the fixed-order expansion of our resummed result from the full fixed-order result. To 
obtain a result which is valid both to NNLL and at NLO in fixed-order perturbation theory. 
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we compute 



dqr 



matched 



dqr 



+ 




expanded to NLO 



(28) 



The qt spectrum is known to 0{al) [19, 20, 21], but it turns out that the matching corrections 
are tiny in the peak region, and 0{as) matching is sufficient for our purposes. To obtain 
the fixed-order expansion of the resummed resuh, we evaluate the hard matching coefficient 
|Cy(— M|, in fixed-order perturbation theory, setting /i/j = ^ in relations (Bl) and (B2) 
of Appendix B, and expand the exponent Qf powers of a^. After this expansion, the 
integrals over the Bessel function in (10) can be computed analytically. A convenient way to 
do this is to first keep the term i]Li^ term in the exponent, and then use (18) to obtain the 
Fourier integral. The Fourier transform of the higher-order logarithmic terms can be obtained 
by taking derivatives with respect to rj of relation (18). The resulting expression was given 
explicitly in equation (60) of [5]. The result for the 0{al) expansion of our result is given in 
Appendix D. 

4 Intercept and behavior near = 

It is an interesting exercise to use the general results of the previous section to derive an 
analytic expression for the intercept of the differential cross section d'^ a / dq^ dy at Qt = 0. On 
one side, this will give us a better understanding of the structure of resummed perturbation 
theory for very small qt- Indeed, we will obtain an expression for the intercept that exhibits 
an essential singularity at = 0, meaning that not even a divergent perturbative expansion 
can be constructed. On the other hand, precise control over the value of the intercept might 
also be of phenomenological relevance in view of future precision measurements at the LHC, 
as it might help to constrain and extrapolate the data into the difficult region of very small 
transverse momentum. It is a curious fact that more than 30 years after the discovery that 
the intercept is calculable provided that the Drell-Yan mass is very large [2], there is still 
no explicit expression in the literature for the 0{1) normalization factor of the result. The 
saddle-point evaluations performed in [2, 4, 22] provide the leading asymptotic behavior for 
Mz oo only, without fixing the normalization of the cross section. Our effective field-theory 
approach allows us to not only derive a closed expression for the intercept at leading order 
in resummed perturbation theory, but we are also in a position to show that corrections to 
this result can be calculated in a power series in ctg. The first-order correction will be given 
explicitly below. 

For Qt = the relevant integral over transverse displacement takes the form (20), where in 
the exponent we must substitute the expression (23) for gp. Using that l — r] = 0{as) = 0{e) 
in this case, we can rewrite this expression as 



9f{v ~ 1, 





a 



'S 



{^0 + A 





:l/2 
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a, d2 + al {kT^^ + 270^^0 + A + 2i^^o + (r^ + A) In ^ ^ 



(29) 



Here ln(g^//i^) counts as 0{1) and allows for scale choices in the vicinity of the default value 
II = q^. The formally super- leading 0{L±) term is absorbed by the change of variables from 
xt^o I — L^, while the leading 0{L\) term needs to be kept in the exponent and provides the 
gaussian weight factor. Expanding out the higher-order terms in (24) and (29) then gives rise 
to gaussian integrals, which can readily be evaluated. In this way, we obtain at next-to-leading 
order in 



Cqq^ij (^1 , , , M| , 11) 



-27E 



27r 



+ as 



as (n + Po) 



1 + Cl{ll) Gs 5{1 - Zl) 6{1 - Z2) Sgi Sgj 



Kl-Z2)5, 



+ 



8(r^ + /3o) 



+ (g,«,^i ^ g,i,^2) 



(30) 



where for Uf — A 



{K + 2-g)^ 
2 (r^ + /3o) 



2 (ro^ + /3o) 



i2(r^ + /3o) 



d2 



^0 + ^ ^2 ctl 
2 ij? 



2 2 
1.48854 + 4.18595 In % + 6.83333 In^ % , 



(31) 



C2(/^) 



2(r,^ + /jo) + 2 ^ V ^ + J' ■ 

Note the peculiar, power-like dependence of the result (30) on the factorization scale. Using 
that for // in the vicinity of we can rewrite 



V^Uas) In 



1 + Xa, + (r^ + /3o) a. In ^ + 0{al) 



cuspv 



(32) 



in (5), it is straightforward to check that the product |Cv(— M|, /i)p Cqq(^ij{zi, z^-, 0, M|, /x) is 
indeed RG invariant. Numerically, the NLO correction proportional to Os in (30) is of rather 
modest size. 

With the default scale choice /j — q*, the prefactor in (30) can be written as 



27r 



/i2 Va,(r^ + 



27re-2T^^ 



TT 



C'fq;s(?*) 



(r^ + ^o) «.(g*)' 



(33) 
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which features an essential singularity at q;^ = 0. When combined with the NLO expression 
for the hard function |Cy(— M|, /x) |^ given in Appendix B, the result (30) provides an explicit 
expression for the intercept of the Drell-Yan spectrum (Pa/dq^dy at qt = 0, to NLO in 
RG-improved perturbation theory. While the leading term is of a genuinely non-perturbative 
nature, higher-order corrections can be calculated in a systematic way in powers of as{Mz) 
(for |Cyp) and as{q*) (for Cqq^ij). Even though the the calculability of the intercept (for 
sufficiently large Drell-Yan mass) is known since the paper [2] , to the best of our knowledge this 
is the first time that explicit expressions for the normalization and the first-order perturbative 
correction have been derived. 

Since resummed perturbation theory allows one to predict the intercept of the transverse- 
momentum spectrum at Qt = 0, it is natural to ask whether similar methods can be employed 
to construct a series expansion of the differential cross section d^a/dq^dy in powers of which 
is valid for qt < q*. To this end, one would like to compute derivatives of the distribution (10) 
with respect to q^ evaluated at qr — 0- Any attempt to do so leads to an encounter with a 
violently divergent series. Taking n derivatives of the Bessel function jQ^XTqr) with respect 
to q^ and setting = generates a factor ~ e"^ in the integrand in (20), and evaluating 
the resulting gaussian integral for ~ 1 yields an extra factor of e" /^"^"-"^^o ^jj^g leading 
term in (30) is then multiplied by the series 

27rn2 

(r^ + /3o) a.(g*) 

(34) 

1 - 1.435^ + 10.66 1^1 - 729.7(^1 +5.82-10M^) 

The resulting divergence is much worse than a factorial growth in n and renders any attempt 
to derive an analytic expression for the shape of the qT spectrum near the origin hopeless. 
Having already seen that the integral over the Bessel function in the first fine of (26) does not 
want to be expanded in powers of a^, we now observe that, even more so, it does not want to 
be expanded in powers of qx- The shape of the distribution for small transverse momenta is 
genuinely non- analytic. 

The discussion just presented raises an interesting issue. In (4), we have presented the first 
term in an OPE of the transverse-position dependent PDFs Sj/at in terms of ordinary PDFs 
defined in terms of the nucleon matrix element of a bilocal quark operator at zero transverse 
separation. In SCET, it would be straightforward to extend this expansion to higher orders in 
^QCD^T- f"*^^ very small qrp we would expect the dynamically generated scale q^, to determine 
the size of these power corrections. However, at any finite order such an OPE will run into 
the same problem as mentioned above, generating "power corrections" whose coefficients grow 
hke e*" In order to avoid this problem the OPE must be resummed, which of course is 
not feasible in practice. On the other hand, reasonable models for such a resummed OPE 
can be obtained by noting that the transverse-position dependent PDFs must vanish rapidly 
when the two quark fields are separated by a transverse distance xt larger than the proton 
size, since then the quark bilinear has very little overlap with the external nucleon state. It is 
thus reasonable to use an ansatz of the form 

Bq/N{i, xl, n) = /hadr(a;rANp) B'^mi^^ /^) > (35) 



V^^f^Vexp 
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where the perturbative functions Bf^^ carry all the scale dependence and are given by (4), 
whereas the hadronic form factor /hadrl?") with /hadr(O) = 1 describes the fall-off at large 
transverse distances and is parameterized in terms of a hadronic scale Anp- For simphc- 
ity, we assume that this form factor is independent of ^. The above ansatz inserts a factor 
[/hadr(3:TANp)]^ uudcr the integral over xt in (10), which suppresses the region of very large 
xt values. We will adopt the models 

/har(a^TANp) = exp (-A^p 4) , /£Sr(a;TANp) = A 2 (36) 

for the form factor, which agree in their first-order terms but have quite different behavior for 
large separation. Fortunately, we will find that while the results are rather sensitive to the 
value of the hadronic scale Anp, the precise shape of the form factor appears to be of minor 
importance. We will show in Section 5 that hadronic corrections to the value of the intercept 
of the da/dq^ distribution indeed scale (approximately) as a power law, ~ {-^qcd/i*)^, where 
due to the resummation of the OPE the exponent 5 is not given by an even integer. Let us 
note for completeness that the hadronic form factor /hadr can in general also depend on the 

A 2 2 

quark flavor and the momentum fraction ^. Models which include a factor ^^2^t in /hadr, with 
a second non-perturbative parameter A2, were studied in [23, 24]. 

5 Systematic studies 

Having discussed the structure of the theoretical prediction for the resummed Drell-Yan cross 
section in detail, we now proceed to perform some systematic studies related to scale variations, 
various implementations of the expansion, different ways to set the scale n, and the importance 
of power corrections. For the purposes of the discussion in this section, we will consider the 
resummed cross section da/dqx without matching to fixed-order perturbation theory. We 
show results at LO and NLO in RG-improved perturbation theory, which correspond to NLL 
and NNLL accuracy. Since we will later match to fixed-order calculations, we will refer in the 
following to the resummed results by their logarithmic accuracy rather than their order in RG- 
improved perturbation theory so as to avoid confusion. For concreteness, we consider the case 
of Z-boson production at the Tevatron with the Z-boson decaying leptonically, — )■ X + Z — )■ 
X + . The corresponding cross section is obtained by multiplying the Z-production cross 
section with the leptonic branching ratio Br(Z — > i'^i") — 0.03366. Throughout, we use 
sin^ 9w — 0.2312 for the weak mixing angle and a{Mz) = 1/128.89 for the fine-structure 
constant. We use MSTW2008NNLO [25] as our default PDF set, which has an associated 
value of as{Mz) = 0.11707. The strong coupling is evolved with three-loop accuracy and has 
fiavor thresholds at iib = 4.75 GeV and /^c = 1.4 GeV for the b and c quarks. 

In Figure 2, we show the scale dependence of the cross section obtained in the improved 
resummation scheme developed in Section 3, varying the factorization scale /x by a factor 
two about the default choice = g^, + qt, where ^ 1.88 GeV has been defined in (19) 
and emerges dynamically in the region of very small transverse momentum. We observe 
a significant reduction of the scale uncertainty when going from NLL to NNLL order. In 
addition to the scale /i, the cross section also depends on the hard matching scale /ih, which 
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Figure 2: Comparison of NLL (blue bands) and NNLL (green bands) predictions for the cross 
section in the improved expansion scheme. The factorization scale fi is varied by a factor two 
about its default value fi = qx + q*, while the hard matching scale is fixed at = — M|. The 
thick lines in the left plot are obtained for the default scale choice. 



2 2 2 2 

NNLL 1.087l°-°J? 1.13llS-[i?l 

N3LL i.imtS 1-1301°:^ 



Table 1: The hard function \Cy{—M'^,fi)\'^ at /i = Mz for space-like and time-like choices of 
fj.\. The uncertainties are obtained by varying fi^ by a factor two about the default value. 

arises when the hard function \Cv{—M'^,fi)\'^ is evolved from a high scale fih ~ Mz to the 
scale /i. The solution of the corresponding RG equation is well known and is reproduced in 
Appendix B. If /Xc < /i < /ib, the hard function has to be evolved across a flavor threshold. In 
this case, we first evolve from fih down to the threshold ^b, switch to four flavors, and then 
evolve from to /x. While the all-order solution is independent of the matching scale fit, a 
residual scale dependence remains at finite orders. To estimate the associated uncertainty, we 
should vary the matching scale in the hard function. However, since the hard function is an 
overall factor multiplying the cross section, we can discuss the fih dependence independently 
from the rest. To separate off the gr dependence arising from the choice of fi, we set fi = Mz 
for the discussion of the /x/i dependence. The perturbative expansion of the hard function and 
its fih dependence are displayed in Table 1, which shows two different choices for the default 
matching scale: the space-like choice fil = M| and the time-like choice fil = — M|. Picking 
n\ = —Mz is motivated by the fact that the vector form factor is evaluated at a time-like 
momentum transfer. In [16], it was shown that this choice greatly improves the convergence 
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qj [GeV] qr [GeV] 

Figure 3: Comparison of different ways to perform the perturbative expansion and scale 
setting. The sohd hne in the two plots shows our default prediction, the dashed curve in the 
left plot is obtained by performing the e expansion in the exponent, while the dotted curve in 
the right plot is obtained if the scale setting is performed for the integrated cross section. 

of the expansion of the Higgs-boson production cross section. In the Drell-Yan case the effect 
is less pronounced [16, 17, 18]. At NNLL order (NLO in RG-improved perturbation theory), 
which is the relevant approximation for our paper, the corrections are similar in both schemes, 
but at N^LL (corresponding to NNLO) the corrections and the associated scale uncertainty 
become much smaller for /i^ = — M|. Using the values in the table, one can easily read off the 
uncertainties associated with the variation of and switch between the two schemes. Our 
NNLL plots are given for /i^ = — M|, so that the hard scale uncertainty is (l?;o)%- To obtain 
the result for /i^ = +M|, one multiplies by 0.961 and obtains a scale uncertainty of 
In view of the known N^LL result for the hard function, it is likely that the value obtained 
with /i^ = — M| is more reliable and that the scale uncertainties obtained with /i^ = M| 
underestimate the higher-order corrections. 

In Figure 3, we compare different ways to perform the perturbative expansion in RG- 
improved perturbation theory. For example, instead of expanding out the higher-order terms 
in e contained m. gp and in the hard function, we can keep these terms in the exponent. The 
result is shown by the dashed line in the left plot. It is gratifying to see that the difference 
between expanding out and keeping the terms is indeed very small in the peak region, such that 
the lines obtained with the two prescriptions are nearly indistinguishable. The small difference 
also justifies treating the hard function as an overall factor. In the numerical implementation 
of the traditional CSS formalism, the scale setting is often performed for the integrated rate 
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Figure 4: Long-distance effects on the differential cross sections (upper row), the intercept of 
da/dq^ (lower left), and the peak of da/dqr (lower right). The solid blue lines are obtained 
from modeling with a gaussian, while the dashed red lines correspond to a pole form. 



and the spectrum is then obtained as the derivative of S(gT)- As long as the renormalization 
scales are held fixed, integrating and differentiating commute. However, if the scale is chosen 
in a gr-dependent way, /i = q^ + qr, there is higher-order scale dependence in which the results 
obtained from I^iqr) differ from a direct evaluation of the spectrum. An advantage of setting 
the scale in E(q'T) is that one automatically recovers the fixed-order result if one integrates 
the resummed spectrum to high values of qx, since the entire integral is then evaluated at 
a high scale. One can argue whether this is an important requirement, since after all the 
resummation was designed for the low qT region, not for the total cross section. Either way, 
the right plot in the figure makes it clear that the associated difference is very small. At larger 
values qr ^ 15 GeV the differences with respect to our default choice, shown in the panels 
below the plots, are no longer negligible. To reduce the scheme dependence in this region one 
will need to match to 0{al) fixed-order results, since the bulk of the difference arises from 
0{a'^) terms beyond the accuracy of our computation. 

Before comparing to data, we now study the effect of long-distance corrections. As dis- 
cussed near the end of Section 4, because of the violent divergence of the expansion of the 
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cross section around = it is not possible to perform the usual twist expansion to study 
these corrections. Instead, we will model long-distance corrections using form factors such as 
the ones shown in (36) under the Fourier integral in (2), where A^p is a typical scale associated 
with low-energy QCD. These form factors suppress the region of large xt- The effect of the 
non-perturbative corrections on the spectrum is shown in Figure 4 for < Anp < 1 GeV. 
Already for > 3 GeV, the effects are almost negligible in da/dq^. It is also remarkable that 
they are to a large extent insensitive to the shape of the form factor. The dashed red lines ob- 
tained with the pole form lie very close to the solid blue lines obtained with the gaussian form. 
For da/dqT, the corrections result in a small shift of the distribution. The dependence of this 
shift on the parameter Anp is shown in the fourth plot in Figure 4 and is quite similar for the 
two different cut-offs. The kink in the peak position arises because the peak happens to move 
over the fe-quark flavor threshold as An? is increased beyond 500 MeV. We integrate out the b 
quark at a scale fih = 4.75 GeV, which corresponds to a qt value of = A^t — 9* ~ 2.9 GeV. 

We have mentioned in Section 4 that we expect the dynamically generated scale to set the 
size of hadronic power corrections to the short-distance results obtained in this paper. With 
the example of the dependence of the intercept of the da/dq^ distribution on the hadronic 
parameter A^p we have tested this assertion numerically, studying a large range of values of 
Anp and the boson mass My (and hence q^,). We have confirmed that the resulting relative 
shifts can be well approximated by a power law ~ A^p/(gf Aq^^), with 6 ^ 1.5 for the 
Tevatron and 5 ^ 1.0 for the LHC with ^/s = 14TeV. 

6 Comparison to experimental data 

We now compare our results to the available experimental data. The most detailed picture of 
the low-gr region is provided by the results by CDF [26] and D0 [27] obtained during Run I 
of the Tevatron, which are quite finely binned at small transverse momentum. In Figure 5, 
we plot the the experimental data, together with our prediction, obtained using the improved 
expansion scheme at NLO, matched to the 0{as) fixed-order result. The small matching 



20 




corrections are shown by dashed hnes in the figure and have been muhiphed by a factor 5 to 
make them visible on the scale of the plots. For the hard scale we choose the time-like value 
/i^ = — M|, and we use /i = g^ + g* for the factorization scale, which is varied by a factor of two 
to obtain the uncertainty bands in the various plots. We do not show the small uncertainty of 
-1%*^° associated with varying the hard scale. It is gr-independent and can be found in Table 
1. The PDF uncertainties are shown in Figure 6 and are of similar size at the Tevatron and 
the LHC At low gr, they are of order 2.5% and then decrease linearly to roughly 1% near 
qt = 30 GeV. In the same figure, we also show the predictions obtained using NNPDF 2.1 
PDF sets [28] . They give somewhat larger uncertainties. The central value is almost identical 
to MSTW at the Tevatron, but about 7% lower in the peak region for the LHC. 

The overall agreement of the data [26, 27] with our result is very good, but the D0 data is 
consistently lower than the CDF result and our theoretical prediction. Summing the data bins 
to the total cross section pp Z + X + X, one finds that the result of CDF amounts 

to (Ttot = 247.4 pb, while D0 obtains Utot = 221.3 pb. The theoretical prediction for the total 
cross section at 0{a^) is atot = 229.7 pb, in between the two values.^ On the right hand side, 
we show the result obtained by normalizing the data to their respective total cross section. 
One observes that the shape agrees well between the two experiments. We have normalized 
the theoretical curve to the 0{a1) fixed order result. Note that this does not guarantee that 
the theoretical prediction for the spectrum is normalized to one, since we have only matched 
to 0{as) in the spectrum. 

In the plots in Figure 5 the observed peak is slightly (by about 750 MeV) to the right of the 
prediction, which was obtained setting the non-perturbative parameter Anp = to zero. We 

^We have used the code VRAP [29], based on the paper [30], to obtain the total cross section. 
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Figure 7: Comparison with Tevatron Run I data from CDF, with and without long-distance 
corrections. The lower panels show the deviation from the default theoretical prediction. 

have discussed in the previous section that long-distance corrections will shift the peak to the 
right, and Figure 4 shows that a shift of 0.75 GeV corresponds to a value of Anp = 0.6 GeV. 
In Figure 7, we compare again to the CDF data [26] and plot the theoretical prediction for 
both Anp = and An? = 0.6 GeV. In the lower panels, we give the ratio of the experimental 
and theoretical results to our default prediction. Including a non-perturbative shift, a good 
description of the data is achieved over the entire qt range. In Figure 8, we repeat the same 
comparison for the Tevatron Run II results from D0 [31, 32] and for the LHC result of the 
ATLAS collaboration [33] . Since this data is not finely binned in the peak region, it difficult to 
draw firm conclusions on the necessity for long-distance corrections. However, in both cases, 
the first data bin is below the prediction without including a long-distance correction. 

The systematic experimental uncertainties which affect the low qt experimental results are 
substantial, because it is highly sensitive to lepton transverse momentum resolution. Recently, 
two new variables ar and 0* were introduced, which probe the same physics but have reduced 
sensitivity to the momentum resolution [34, 35]. D0 has now performed a very precise mea- 
surement of the variable 0* [36]. It would be interesting to include the lepton decay in our 
results and to study these variables. In the traditional framework, resummed results for these 
quantities were presented recently in [37, 38]. 

The region of larger qx > 20 GeV is not affected by long-distance corrections and should be 
described well by fixed-order perturbation theory. In this region the data lies somewhat above 
the prediction, in particular for the case of the ATLAS results. A comparison to the existing 
fixed-order results is given in Figure 9. The red bands correspond to the 0{a1) fixed-order 
result for the spectrum, which the highest order currently known. To compute this result we 
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Figure 8: Comparison to Tevatron Run II and ATLAS data, with and without long-distance 
corrections. The lower panels show the deviation from the default theoretical prediction. 



use the numerical code QT [39]. For the sake of comparison, we have evaluated all results 
using the MRST2008NNLO PDF set. The fixed-order results diverge to ±00 for vanishing 
transverse momentum. Since the fixed-order result depends both on qt and Mz it is not clear 
which value one should choose for the renormalization and factorization scales. The edges of 
the fixed-order band in Figure 9 correspond to the two choices fi = qt and /i = Mz- The 
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Figure 9: Comparison of the resummed result (green) and 0{a1) fixed-order result (red) for 
qt spectrum at the Tevatron and at LHC 



choice n = qt diverges to — oo, while = Mz rises to +00. Our result is consistent with 
the 0{a1) fixed-order result at large qx, but performing the matching to this order would 
presumably lead to a slight increase of our cross-section prediction in this region. 

We have also compared our findings to the results [40]. They perform resummation at the 
same level of accuracy but with a different formalism, based on the CSS formula [4]. The 
peak position is exactly the same as in our result, and their peak height is about 2.5% lower, 
but compatible within uncertainties. Since we have adopted a time-like hard matching scale, 
which gives a 4% higher cross-section than the usual space-like choice, this difference is not 
unexpected. Also at qt = 30 GeV, the two results are compatible, but their central value is 
about 8% higher than ours, perhaps because the NLO matching has been implemented, or 
because of the unitarization prescription implemented in their formula. By inspecting plots 
showing a comparison of results of the RESBOS code [41] with Tevatron [42] and LHC [33] 
data, we conclude that our results are compatible with these predictions, which also have 
N^LL accuracy. We also note that our results have much smaller scale uncertainties than the 
numerical results obtained in [14, 43] based on the formalism presented in [13]. Figure 3 in [43] 
shows that the peak height in da/dq^ changes roughly by a factor 3 when the factorization scale 
II is varied in the narrow range 2.8 GeV < /z < 4.4 GeV. In our case, the peak height changes 
only by 5%, even though the scale is varied over the larger range 2.1 GeV < fi < S.6 GeV. We 
believe that this numerical instability illustrates our earlier point, namely that the results of 
[13, 14, 43] only have LL accuracy (in the exponent), since they do not resum the corrections 
associated with the coUinear anomaly. 
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7 Conclusions 



The transverse momentum spectrum in Drell-Yan processes is one of the most basic observ- 
ables at hadron colUders. It neverthess manifests a number of remarkable properties at low 
transverse momentum. These are associated with the collinear anomaly, which complicates 
the factorization of the cross section. Because of the anomaly, not only the hard function, 
which encodes the virtual effects associated with the electroweak boson production, but also 
the product of collinear functions has a dependence on the boson mass My. In a previous 
paper, we have shown that the anomalous My dependence exponentiates in transverse po- 
sition space. This guarantees that the associated large logarithms are under control, but at 
the same time renders the Fourier integral non-trivial, which yields the transverse momentum 
spectrum. One of the interesting consequences of the presence of the anomaly is that the scale 
which controls the perturbative corrections to the factorization theorem is no longer /i qt 
at very low qt, but saturates to a non-perturbative value q* ~ Myc''^"^^^^"'^^^^ for qt 0. 
Numerically, q* ~ 1-88 GeV for Z-production, which implies that the spectrum remains short- 
distance dominated even in the limit qT 0. The fact, that the zero intercept of da/dq^ 
is calculable for very large Drell-Yan masses was pointed out a long time ago by by Parisi 
and Petronzio, but our formalism has allowed us to systematically evaluate the corrections 
to their asymptotic result. The uncertainty on the calculable short-distance corrections, as 
well as the long-distance effects, which we estimate using models, are seizable given that the 
relevant scale is q* ~ 1.88 GeV. Nevertheless, it is possible to obtain a prediction for the in- 
tercept for Z-production with an accuracy of about 25%. This demonstrates that, contrary to 
common belief, the mechanism identified by Parisi and Petronzio is not merely an interesting 
theoretical result, but of phenomenological relevance. 

In addition to making the spectrum infrared safe, the anomaly precludes the expansion 
of the spectrum in any of its parameters. In particular, there is a strong factorial divergence 
in the expansion in ag, associated with terms in the Fourier integral which become enhanced 
near g-^ = 0. While they are not logarithmically enhanced, these terms nevertheless need to 
be resummed. An even more violent divergence arises when one tries to expand the spectrum 
around g|, = 0. In this case the n-th order terms in the expansion are enhanced by e"" , which 
renders the expanded result meaningless. To study long-distance effects, we insert a cut-off 
function such as a Gaussian exp(— A^pxf.) into the Fourier integral. We find that the long- 
distance effects are of the expected size and largely independent of the form of the cut-off. 
However, a twist expansion exp(— A^pXy) = 1 — A^pxf. + . . . suffers from the same strong 
divergence as the expanision around g^- = and can thus not be used. 

We have compared to experimental results from both runs of the Tevatron and from the 
LHC. We find excellent agreement with our predictions, however, for ^ 3 GeV the agree- 
ment is only achieved after including long-distance effects. The same choice of the associated 
parameter Anp = 0.6 GeV describes these effects both at the Tevatron and the LHC, which 
provides a consistency check on the model we used. It would be interesting to study the long- 
distance effects in more detail. To this end, the recently proposed variables Gt and 0* would 
be especially well suited, and the variable (p* has now been measured by DZero. To compute 
these quantities, one needs to include the decay of the electroweak boson into the lepton pair. 
It is straightforward to extend our results to this case, and it is important in order to be able 
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to implement the lepton cuts used by the experiments, and to compute the charged-lepton 
spectrum used for the W^-mass determination. 

The formahsm developed here can be applied to other Drell-Yan-type processes, such as 
Higgs-boson production or the production of new heavy, color-neutral particles at hadron 
colliders. The case of Higgs production via gluon-gluon fusion is particularly interesting. In 
this case many of our formulas apply with a simple substitution of color factors. In particular, 
in the expressions for the dynamically generated scale in (19) and for the exponent gp of 
the Bessel integral in (23), one must replace the one-loop cusp anomalous dimension in the 
fudamental representation with that in the adjoint representation, Fq — )■ Fq . As a result, for 
a Higgs with mass of 120 GeV, for example, the value of is 7.5 GeV. This is truly a short- 
distance scale, and all of the methods developed here should work much more accurately than 
in the case of Z-boson production. In particular, long-distance effects in the region of very 
small qT should be strongly suppressed. A phenomenological analysis of Higgs production in 
our framework is left for future work. 
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A Single differential cross section 



Integrating the double differential cross section (1) over rapidity, we obtain 



da An'^a 



dqt 



gv(-Mi,/.)rE ''^''t''''\ ^s:/ 



'^dz 
z 



i=<l,9 3=1,9 



(Al) 



X 



Cqq^ij {z, q^, M|, //) ffij{r/z, ij) + (q,i^ q, j) 



where g^ = (g/e) {T^ — Cg sin^ 6\v) and gf^ = (g/e) (— Cg sin^ 9w) denote the couplings (in units 
of e) of the Z boson to left-handed and right-handed quarks, and hence 



biP + \9r\^ _ (1 - 2|e,| sin^ Ow)^ + sin^ 



8 sin cos^ Qw 



(A2) 



The parton luminosities are defined as 

du 



u 



^i/m (U, ll) (pj/N^iz/u, ll) = (^(pi/Ni ilj) ® (pj/N2 ()")) (z) ■ (A3) 
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The kernel functions Cgg^ij are obtained from tlie functions Cgg^ij defined in (2) by a convo- 
lution in the two Zi variables, such that 



1 r / ^2 j\ .f2\ -P'qqiL±,a,a) 



-'0 



X 



,as)<S) Iq^j{L±,asj^{z) . 



(A4) 



B Hard matching coefficient 



The solution to the RG equation (5) for the hard function takes the general form [15] 

(Aj2 \ -2ar(M?i,/i) 
^ |Cy(-M|,/i,)|', (Bl) 

f^h J 

where ~ is a hard matching scale, at which the value of Cy is calculated using 

fixed-order perturbation theory. At one-loop order 



|Cy(-M|,/.,)f = 1 + 



27r 



Re -L^ + 3L - 8 + 



+ 



(B2) 



where L = ln(— M|//x^). The two-loop correction can also be found in [15]. The advantages 
of using a time-like scale choice {fi\ < 0) for time-like processes such as Drcll-Yan production 
were emphasized in [16, 17]. The Sudakov exponent S and the exponents a„ are given by [44] 



S{l^h,IJ') = 



^ 



An 



1 1 1 ^^^rf A 



- I (1 — r -|- Inr) H In^ r 



+ 



An 



rf/3i 132 

ro^/3o /So 



[1 — r + rlnr) + 



Pi ^2 



PI P2 rf/3i 



+ ^ 



PI Po 
1 -r)2 



(1 — r) Inr 



(B3) 



^ 



Pi Po r^Po 



+ . 



where r — as{(J^) / as{(^h) . A similar expression, with the coefficients Ff replaced by jf, holds 
for the function a^q. The relevant expansion coefficients of the anomalous dimensions and 
^-function can be found, e.g., in [15]. 

For the special scale choices fih = Mz and iJ = q^, the expression for the hard matching 
coefficient can be simplified, since the two couplings as{Mz) and as{q*) are related via the 
condition that rj — 1. Solving the relation 



as{Mz) 



as{q*) 



da 

W) 



1 , M| 
- In^ 

2 ql 



2F^a, 



(B4) 
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with as = Q;s(g*)/(47r) iteratively, we derive 



ocs{Mz) 



1 + c 



+ a, 



/3o 



In 



1 + c 



+ 



1 + c 



Pi 



cln 



1 + c 



+ 0[al{^i)], (B5) 



where we have introduced c = / = 16/25. For the purposes of our discussion in this 
appendix we neglect the presence of flavor thresholds and use Uf — A active quark flavors 
throughout. Using the above relation, we find that 



|C^(-M|,g,)| =exp 
where 
jVf = exp 



1 - (1 + c) In 



1 + c 



1 + 1.51591 a^ + C>(a^) , (B6) 



A 



cln^ 



1 + c 



1 

Jo 



rf 



1 - cln 



1 + c 



= 1.83152. 

(B7) 

The explicit expression for the 0{as) terms is more lengthy, and we refrain from giving it here. 

Combining the above result for the hard matching coefficient with expression (30) for the 
coUinear kernel functions at qt — 0, we obtain for the default choice /i — q* 



-27B 



exp 



Potts 



+ f lln 



1+* 

^ ^ -pF 
^ 



27r 



{T^ + Po)as 



(B8) 



X ^Vf { (1 + 3.00445 as) S{1 - zi) S{1 - z^) Sgi % + ...}, 



where the remaining terms in the parenthesis have the same form as in (30), evaluated at 
/i = g*. Note that the non-perturbativc exponential factor in the first line can be expressed as 
e^o^ L - J _ j ■ Along with (1), formula (B8) provides an explicit expression for the 

intercept of the Drell-Yan transverse- momentum spectrum in terms of the scale g*. By using 
the one-loop expression for the running coupling as{ij) to eliminate in favor of the ratio 
M^/Aqcdj we can relate our result to approximate saddle-point expressions for the intercept 
derived in [2, 4, 22]. We find that (B8) is consistent with relation (58) in [22] up to subleading 
terms involving the two-loop /3-function. However, the correct normalization factor Mf is 
derived here for the first time. 



C Perturbation theory with modified power counting 

For a consistent evaluation of the cross section including 0{e) corrections in our modified 

power- counting scheme, where Lj^ ~ 1/ ^Ja~s, we need to retain some terms of up to four-loop 
order in F^^, up to three- loop order m Hf, and up to two-loop order in the kernels Iq^i. The 
relevant terms are, however, determined by the RG in terms of known anomalous dimensions 
and ^-function coefficients. We obtain them by solving the evolution equations (6) and (9) 
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recursively, using that d/d\n/j = 2d/dL± + P{as) d/das [5]. The terms that need to be 
retained for the purposes of our work read 



Fqq{L±, ttg) — tts Tq L±_ + 



r^^o^ ^ + (ro^/3i + 2rf /3o) ^ 



Po ^ ^ ■ ■ ■ 



+ ... , 



hF{Lj_,as) = as 



r^/3o 3^ + (rf - 27oVo) ^ + . . . 



i OPO + • • • 
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+ ..., 



(CI) 



75^1(2;, L_L, a^) = 5(1 - z) dqi + a 

(Vqaz)-2/3oVl'l,{z))^ + ... 



+ 



The expressions for Iq^i involve the convolutions Dq^i of DGLAP splitting functions defined 
in (25). Using the one-loop splitting functions in (16) along with 



Vflq{z) = 8C, 
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we obtain 
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In deriving the first result, we have used that 

1 1 



(i-z)+ {i-zY 



In 



(C4) 



D Matching to fixed-order perturbation theory 

Here we evaluate our resummed expressions for the hard matching coefficient \Cy\ and the 
kernels Cgg^ij in fixed-order perturbation theory, in terms of the coupling = Q;s(A'')/(4vr). 
We choose jj, ~ qt, so that Lm = ln(M|/g|,) is the large logarithm, whereas — ln(/x^/g|^) 
counts as an 0{1) number. 

The one-loop, fixed-order expression for the hard-matching coefficient can be taken from 
(B2). We find 



|Cv(-M|,/.)f = l + a, 



7-^ 

— ^ {Lm - L^Y ~ 270 {Lm - L^,) + Cf[ — 



- 16 



+ 0{al) . 



(Dl) 

There is no need to include higher-order terms in this expression as long as we restrict ourselves 
to the spectrum at qr 7^ 0. To obtain the fixed-order expansion of the kernel functions, we 
start from our resummed expressions obtained using the e expansion and reexpand them to 
second order in a^. This gives 

ql Cgg^,,j{zi, Z2, ql, Ml, ji) 

as {TqLm + 27o) + al ki{LM, L^) 8{1 - zi) 8{l - Z2) Sgi Sgj 
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^ + al k2{LM, L,)] [vl'U^,) 5(1 - Z2) % + d{l - z^) dg, V^l^{z2) 



+ «S { ^qLm TZg^i{Zl) 5(1 - Z2) 6gj +6{1- Zl) Sgi ng^j{Z2) 



Vq^i{Zi) 6{1 - Z2) Sgj + 5(1 - Zl) 5giVq^j{z2) 



(D2) 



where 

ki{LM, Ln) 

k2{LM, Lfj) — 



L, {T^Lm + 27o') 



+ rfLM + 2(r^)\3, 



-L, (r^LM + 27o^) + ^lLj + ^L,, 



(D3) 



We emphasize the important fact that, due to the coUinear factorization anomaly discussed in 
Section 2, the kernel functions Cgg^ij contain a dependence on the hard scale Mz via the large 
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logarithms Lm = 1ii(^I/Q't)- These large logarithms are not resummed through the evolution 
of the hard matching coefficient |Cv|^, but as shown in (2) and (10) they exponentiate in xt 
space. At n*'* order in perturbation theory, the Born-level structure 5qi Sqj in (D2) receives 
corrections of order {agLM)"', which are of 0{1) in RG counting and must be resummed to 
all orders to accomplish a consistent resummation with NLL accuracy. As we have already 
mentioned in Section 2, the approach of [13] fails to resum the large logarithms in the kernel 
functions Cqq^ij (which are written in the form Q = Xn®Tn® in this paper). The one-loop 
logarithmic term TqLm in the first structure in (D2) was reproduced in this approach and 
it resulted from the convolution of coUinear beam functions with an inverse soft function. If 
it could be shown that, to all orders in perturbation theory, only a single large logarithms 
arises from SCET convolution integrals, then indeed there were no need to resum them and 
the approach of [13] would be consistent. However, already the presence of the two-loop 
logarithmic term — (r^)^L|^L^ entering via the coefficient ki in (D3) shows that this cannot 
be the case. The presence of this term was missed in [14] by making the scale choice ji — qt, 
which sets = 0. With any other scale choice a quadratic term ~ (ilL\j is required by RG 
invariance. The peculiar feature that this term vanishes for ^ = q^is due to the fact that in the 
MS scheme the anomalous exponent for Drcll-Yan production does not contain a non-trivial 
constant at one-loop order, Fqq{L±, as) = as {TqL± + dl) with df = 0, sec (CI). In any other 
scheme except the MS scheme, the coefficient ki contains a term —dl TqEI^. We also note that 
in other applications of the collinear anomaly, such as in the study of resummation for the jet 
broadening distribution in e'^e~ annihilations [12], the corresponding one- loop coefficient di 
is non-zero. 

Adding the contributions from the hard function, we find that the product \ Cv{—Mz, 
xCqq^ij{zi, Z2, qT, IJ>) obeys the same decomposition as shown in (D2), with coefficients 
^1,2 replaced by ki^2 given by 

A;i(Lm, L^) = (r^ Lm + 27o') - 27o' Lm + PoL^ + Cp f — - 16 

+ rfLM + 2(ro^)'C3, m 

Note that in this result only single logarithms remain. It is a remarkable fact that with 
our improved expansion scheme based on the modified power counting these scale-dependent 
logarithms are precisely those ensuring that the differential cross section is scale invariant 
through 0{a^). This did not have to be the case, since at two-loop order L^ terms not 
accompanied by a large logarithm Lm are in principle beyond the accuracy of our NNLL 
calculation. For example, performing the resummation based on a conventional expansion in 
powers of would miss the last two lines of (D2) as well as some other //-dependent terms. 
The fact that with our modified power counting we do account for all scale-dependent terms 
at 0(0^) explains why in this case the bands are much narrower than those obtained with a 
conventional power counting, see the red and green bands in Figure 1. 
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